Mon. Not. R. Astron. Soc. OOP. [TIIT^ (2005') Printed 5 February 2008 (MN WF^ style file v2.2) 



A model for the metallicity evolution of damped Lyman-a 
systems 



P.H. Johansson* and G. Efstathiouf 

Institute of Astronomy, Madingley Road, Cambridge CBS OHA, UK 



Version of 5 February 2008 



o 
o 

(N 



(N 
> 

cn 
o 

o 
o 



X 



ABSTRACT 

We apply a physically motivated stellar feedback model to analyse the statistical 
properties of damped Lyman-a systems (DLAs) expected in the concordance cold 
dark matter (CDM) model. Our feedback model produces extended low-metallicity 
cold gaseous discs around small galaxies. Since the space density of galaxies with low 
circular speeds is high, these discs dominate the cross-section for the identification of 
DLAs at all redshifts. The combined effects of star formation, outflows and infall in 
our models result in mild evolution of the A^Hi-weighted metallicity content in DLAs 
with redshift, consistent with observations. According to our model, DLAs contribute 
only a small fraction of the volume averaged star formation rate at redshifts z ^ b. Our 
model predicts weak evolution in JIhi over the redshift range z = — 5. Furthermore, 
we show that the cosmological evolution of fini and the cosmic star formation rate are 
largely disconnected and conclude that the evolution of JIhi as a function of redshift is 
more likely to tell us about feedback processes and the evolution of the outer gaseous 
components of small galaxies than about the cosmic history of star formation. 

Key words: methods: analytical - intergalactic medium - quasars: absorption lines, 
- galaxies: formation - galaxies: ISM - cosmology: theory 



1 INTRODUCTION 

Damped Lyman-a systems (DLAs) are defined as quasar 
absorption systems with a column density of A^hi J? 2 x 
10^°cm~^, i.e. with neutral hydrogen column densities sim- 



ilar to present day galactic discs (see e.g. IZwaan et alJ 
l2005lJ) . DLAs probe the high column density end of the 
distribution of absorption line systems and are particularly 
interesting because their high concentrations of neutral hy- 
drogen suggests that there may be a close connec tion with 
disc galaxies that we see today JWolfe et alJIlQSd) . In addi- 
tion observations of DLAs measure, in a model independent 
way, the total neutral hydrogen mass density of the Universe 
as a function of redshift. 

There are two main competing scenarios for the ori- 
gin of the DLAs. In the first model, the DLAs are inter- 
preted as large hi gh-redshift progenitors of present day mas - 
sive spiral discs IWolfe et alj Il986t iLanzetta et alJ Il99ll) . 
These discs are assumed to have formed at high (z ^ 5) red- 
shift and to have evolved little apart fr om converting their 
gas in to stars. In support of this theory IProchaska fc Wol^ 
il998h argue that the kinematics of metal absorption lines 
in DLAs are best explained by models of rapidly rotat- 
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ing {vc Xi 200 kms~^) galactic discs. The alternative pic- 
ture, which is more in line with current views of hierar- 
chical structure formation, posits that the DLA systems 
are 'dwarf galaxies rather than fully formed ^ L* disc s 
jKauffmann fc Charlotill994l: iMo fc Miralda-Escudlll^^4^ . 



In particular. "^a ehnelt et al.l il998D used hydrodynamic 
simulations to show that the statistical properties of the ve- 
locity structure seen in the metal lines could be reproduced 
by infalling sub-galactic clumps within dark matter haloes 
with low virial velocities of (fvir ~ 100 kms^^). 

On the observational side, the nature of the DLA host 
galaxies at high redshift is l argely unknown, w ith only a 
handful of optical detections. M ailer et al.l ll2002ll identified 
three DLA host galaxies at z ~ 2 — 3 and found them to have 
luminosities much fainter than L* . At lower redshifts, z ^ 1, 
approximately 30% of t he DLA host galaxies have luminosi- 
ties larger than L* (see IZwaan et alJl2005br) . However, it is 
dangerous to extrapolate from low redshift to high redshift 
since it is plausible (perhaps likely) that the cross-section 
based selection criteria select different objects at different 
redshifts. 

The first comp rehensive metallicity surveys of DLAs 
IPettini et al.l lll994t) derived a global mean A'ni-weighted 
metall icity of logfZ/Zp)) » — 1 at g ~ 2. Further o bserva- 
tions dPettini et al.lll999t |p7ochaska fc Wolfell2000D found 
similar metallicities and no evidence for evolution of the 
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metallicity in the redshift range z ~ 1 — 3.5. At lower red- 
shift {z ^ 1) the situation is more uncertain because of the 
small sample o f DLA s. Recent measurements at z 1 by 
iKulkarni et al.l (|200^ indicate that the evolution of metal- 
licity for DLAs remains weak at lower redshifts and does 
not rise up to solar, or near-solar, values hy z — 0. In ad- 
dition, the cosmological density of neutral hydrogen (f2Hi) 
as traced by the D LA population s hows little evolution be- 
tween z « 0.5 - 5 iRao et alj|200d) . 

Attempts to model DLAs in cosmological simulations 
have proved difficult because of the high numerical res- 
olution required and the difficulties in modelling stel- 
lar feedback and ga lactic winds(see e.g. Katz ct all Il996| : 
iGardner et alj Il997l . |2001^ . In addition to limited numeri- 
cal resolution, early simulations did not follow the chem- 
ical evolution of the DLAs in any detail. Some of these 
shortcomings have been overcome with recent simulations 
by iCcn ct al . ( 20ol) and iNaaaminc ct al. ( 2004b). Both 
studies are able to reproduce the flat metallicity evolu- 
tion of the DLAs, but their metallicity values tend to lie 
at around log{Z/ZQ) ~ —0.5, considerably lower than ob- 
served. To reconcile this higher value with the observa- 
tions, both groups invoke a large dust obscuration correc- 
tion. A dust bias has also been invoked to explain the 
lack of DLA absorbe r s wit h [Zn/H] + log[N(HI)] > 21 
noted b y ^oisse e^ahl l|l998ll . However, a recent study by 
iMurphv &i Lis ke 7200?) using a large sample of quasar spec- 
tra from the SDSS Data Release 2 did not find any evidence 
for dust-reddening caused by intervening DLAs at z ~ 3. 
In addition DLAs extracted from samples of radio selected 
QSOs, which should be largely free from any dust bias, show 
only marginally higher metallicities co mpared to absorbers 
from optically selected control samples iEUison et alj I2OOII : 
lAkerman et al.ll2005l : Ijorgenson et"al]l200d) . 

Semi-analytic modelling of DLAs provides a comple- 
ment ar y approach to num eric al s im ulations (K auff mannI 
| l996l: iPrantzos fc Boissicr 2 00o|; iMaller et all l200lt 
IsomerviUe et alJ 12001 : Okoshi etalJ l2004) . The semi- 
analytic approach is, of course, much faster than numerical 
simulations and so can be used to explore parametric 
representations of physically complex processes such as 
stellar feedback. Semi- analytic modelling can also provide 
physical insight, which is sometimes difficult to acquire from 
numerical simulations. The disadvantage of semi-analytic 
modelling is that some aspects of the models may be 
sensitive to processes that are extremely difficult to model 
analytically. Even if this is the case (as it is in the models 
described in this paper) it is useful to identify aspects of a 
theoretical model that are robust, and can be understood 
analytically, and those that are likely to require further 
investigation using numerical simulations. 

In this paper we develop a semi-anal ytic model o f DLAs 
using the feedback model described by lEfstathioul (|200(J) 
(hereafter EOO). Some aspects of the application of this feed- 
back model to DLAs were sketched in EOO. The purpose of 
this paper is to develop a more detailed model that can be 
compared with a wide range of new observations. The EOO 
feedback model produces extended cold gaseous discs around 
'dwarf galaxies (defined crudely in this paper as any sys- 
tem in a halo with a circular speed Vc 100 kms~^) . In 
CDM-like models such gaseous discs would dominate the 
cross-section for the identification of DLAs at high redshift 



because the space density of haloes with low circular velocity 
is high jKaufi'mann fc Charlodll994l: IMo fc Miralda-Escud j 
199'!). our model, most of the cross-section is dominated 
by largely unprocessed gas in the outer parts of the galaxies 
that does not participate in the bulk of the star formation 
process. As a result, the evolution of flm in this model is 
linked only indirectly to the star formation history. There 
are some similarities between ou r models and recent semi- 
analytic modelling of DLAs bv Okos hi et alJ ^004). The 
main differences, which are particularly significant for mod- 
elling the metallicities of DLAs, are: (i) our stellar feedback 
and outfiow models are physically motivated, i.e. they are 
related to a well defined, though simplified, model of the 
interstellar medium in galaxies; (ii) rather than applying 
ad-hoc rules of star formation and chemical enrichment, our 
star formation prescription is based on the self-regulated^ 
model described in EOO. 

The structure of the paper is as follows. Section |21 sum- 
marizes the feedback model and its application to DLAs. 
In Section 1^ the feedback model is applied to calculate the 
metallicity evolution of the DLA population. Theoretical 
predictions of various global properties of the DLA popu- 
lation are described in Section^] In Section|^we summarize 
our results and present our conclusions. 



2 THE MODEL 

The EOO feedback model is used to evolve a grid of galaxy 
models as described below. In this Section we will summarize 
briefiy some aspects of the model and describe the relevant 
parameters involved. We refer the reader to EOO for a more 
detailed description of the model. All calculations in this 
paper assume a concordance CDM cosmology with param- 
eters^ 0^ = 0.27, = 0.73, D.t = 0.044, h = 0.71, e.g. 
JSpereel et al. 2003). We assume a scale-invariant adiabatic 
fluctuation spectrum, normalized so that ug — 0.84. 

2.1 Feedback model 

EOO developed a self-regulated feedback model in which the 
instantaneous star formation rate is regulated by disc in- 
stabilities and in which supernovae blastwaves convert some 
fraction of the cold disc gas into a hot component which can 
flow out of the syst em. T he model is based, in part, on the 
iMcKee ~Ostrike r I^T^) theory of the interstellar medium, 
but also includes the infall of gas within a dark halo and the 
formation of a stellar disc through self-regulated star forma- 
tion. 

The dark matter halo of the model ga l axies is described 
by a static NEW profile jNavarro et aljlT99^ . A galaxy 
forms through infall of gas from the surrounding halo un- 
der the assumption that the angular momentum of the disc 
material is approximately conser ved during the collapse of 
the disc jFall fc Efstathioulll98dl) . The angular momentum 
distribution within a dark matter halo is based on a fit to 

^ There is considerable empirical sup port in favour o f self- 
regulated star formation in disc systems ^Kennicutdll998^ . 
^ where ftm, and fij, are the cosmological densities of matter, 
dark energy and baryons respectively and h is defined such that 
Ho=100h kms-^Mpc-i 
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Figure 1. The left-hand panel shows the evolution of hydrogen gas surface density and the right-hand panel the evolution of the gas 
metallicity for a model with A = 0.075, "Uvir = 100 kms~^. The results are plotted at ages of 0.4 Gyr (solid line), 0.8 Gyr (dotted line), 
1.4 Gyr (dashed line), 2.6 Gyr (dash-dotted line), 5.2 Gyr (dash-triple dotted line) and 10 Gyr (long dashed line). The horizontal dotted 
line in the left-hand panel shows the DLA selection limit of 1.6 Mqpc~^. 



the results of A^'-body simulations and is normalised to re- 
produce a target value of the dimensionless spin parameter 
Xh = J\E\^'^G~^M-^''^ , where J, E and M are the angu- 
lar momentum, binding energy, and mass within the virial 
radius ry. 

The stability of the resulting gas eous disc is determined 
by a lGoldreich fc Lvnden- Belli lll965l) like criterion, 

a,^^g{a,p), (1) 
ft 

where jig is the gas surface density, k is the epicyclic fre- 
quency, a and 13 are defined as a* = aag, /x* = Pfig and 
(ji(a, 13) is a correction factor of order unity that applies to a 
two-component gaseous- stellar disc (see EGO for details and 
Fig 1. in EGO for a plot of g{a, 0) as a function of a and P). 

Following .McKee fc Ostrik er (1977J the interstellar 
medium of the galaxy is modelled as a multiphase medium. 
Most of the gas is assumed to be in cold clouds. The system 
of cold clouds is assumed to be marginally unstable and so 
their velocity dispersion is fixed in terms of the surface mass 
densities of gas and stars by equation Q. Supernova blast- 
waves propagating through the interstellar medium convert 
some of the cold clouds into a low-density hot component. 
The instantaneous star formation rate is fixed by balancing 
the energy dissipated in collisions between cold clouds with 
that supplied by superno vae shells . 

Assuming a standard lSalpete'il (Il955l) stellar initial mass 
function (IMF) with mass cutoffs of m\ = O.IM0 and = 
5OM0 and that each star of mass greater than 8Mq releases 
10^^ £51 erg in kinetic energy in a supernova explosion, the 
energy injection rate is related to the star formation rate by 

£sN = 2.5 X lO^'^SsiM* erg s"\ (2) 

where M* is the star formation rate in Mq per year. A 
fraction of this energy rate, ecEsN is assumed to balance 
the rate of energy loss per unit surface area due to cloud 
collisions 



Etii = 5.G X 10''' (1 + ^) ^l%(Jg^ erg s-'pc-^ (3) 

where /igs is the surface mass density of the gas component 
in units of 5AfQpc~'^ and (jgs is the cloud velocity dispersion 
in units of 5 kms~^. 

The uncertain parameter is fixed by normalizing the 
typical star formation rate for a Milky Way type galaxy. 
Assuming a fiat surface mass density profile for the gas to 
J?max = 14 kpc and /3 « 10, a « 5 results in 

ecM* = 0.G04. (4) 

The choice ec = O.Gl results in a reasonable net star forma- 
tion rate of 0.4 MQyr~^ for a Milky Way-type galaxy. We 
therefore adopt = G.Ol throughout this paper. We also as- 
sume a cooling function for gas of a primordial composition 
with a sharp l ower cutoff at T = IG^K. 

Following Mc Kee fc Ostrikeil il977l) an expanding su- 
pernova remnant will evaporate a mass of 

Mev « ■illEl{''{AT^aiN,i(i,^f'''nl*'H4Q (5) 

where rih is the density interior to the supernova remnant 
(in units of cm~^), a; = G.5 pc is the lower limit of the 
distribution of cloud radii and A^d is the number density 
of the clouds (in units of pc"''). (Note that for the typical 
sizes and densities of the cold clouds, the filling factor for HI 
absorption will always be greater than unity for HI column 
densities above the DLA threshold) . The second free param- 
eter of the model, 0k, quantifies the effectiveness of clas- 
sical thermal evaporation (via conductivity) of the clouds 
Kcff = i^i'K S'lid is expected to be less than unity if conduc- 
tivity is reduced by tangled magnetic fields, turbulence, etc. 
The parameter 0k controls the strength of stellar feedback. 
In our normal model we adopt 4>k = G.l, but we have also 
run models with weaker feedback 0k = G.Ol and stronger 
feedback 0k ~ 1.0, (WFB and SFB models, respectively) so 
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that the reader can gauge the sensitivity of our results to 
uncertainties in the strength of stellar feedback. 

The cold gas component is assumed to be in the form of 
HI, which is clearly an oversimplification since some of this 
gas could be converted to molecular hydrogen. However, this 
is likely to occur in the central high-density regions of galax- 
ies rather than the outer parts that dominate the DLA cross- 
section. This together with the low molecular ga s fractions 
observed for DLA systems (sec Wolfe ct al. 2005 and refer- 
ences therein) suggests that ignoring the molecular compo- 
nent will not affect our results. We also ignore the meta- 
galactic UV flux since at the column densities of the DLA 
systems they are of course self-shielded to photoionization. 

The hot phase of the interstellar medium can escape the 
galaxy if the wind speed Ww exceeds the escape velocity Vcsc 
from the centre of the galaxy (see EGO for details on how 
Ww is calculated). A wind with Uw < Wosc is assumed to par- 
ticipate in a galactic fountain and is therefore returned to 
the disc on a dynamical timescale, as described in EGO. The 
model thus incorporates simultaneous inflow and outflow of 
gas ignoring any interactions between these flows. As dis- 
cussed in EOO this may not be a bad approximation since the 
feedback may produce a mildly coUimated outflow with infall 
occurring primarily in the equatorial plane, this type of be- 
haviour has been s een in numerical s i mulations with stellar 
feedb ack (see e. o. iTassis et alj booj : ISoringel fc Hernauisli 
l2003^ 

Chemical evolution is included using the instantaneous 
recycling approximation. Separating between primordial in- 
falling gas accreting at a rate dfii /dt with metallicity Zi and 
processed gas from the galactic fountain of metallicity Zf 
accreted at a rate dfip/dt results in 

HgdZ — pdfis + {Zi — Z)dfj,i + (Zf — Z)df^F- (6) 

Here p = G.02 jHirschi et al.1 l2005ll is the yield and, un- 
less otherwise stated, we assume that primordial gas has 
zero metallicity. The metallicities are normalized to the so- 
lar valuefor which we adopt the new value of Zq = 0.G133 
bv lLodderd ll2003ll . 



2.2 Modelling the DLA population 

Using the model described above we computed the gas and 
stellar surface densities, and their respective metallicity pro- 
flies, for a grid of models parameterized by the virial veloc- 
ity of the parent dark halo (deflned as the circular speed 
at the virial radius) and the spin parameter A. We com- 
puted a grid of 17 models with virial velocities in the range 
Vvii = 35 — 350 km/s and 7 spin parameters in the range 
A = 0.G21 - G.15 resulting in a total of 119 models. The 
lower cutoff of the virial velocity was chosen because gas in 
haloes with low circular speeds will be heated by the ex- 
ternal UV background at high redshifts and so is unable 
to cool (Recs 1986; Efstathiou 1992). (The precise value of 
this lower cutoff is unimportant for most of the results pre- 



sented in this paper. Over the redshift range , 



5, the 



smaller cross-section of such low circular speed systems more 
than compensates for their larger space density. As discussed 
in Section 4, the peak contribution to the DLA population 
occurs at higher circular speeds.) Radial profiles for each 
model were stored for thirteen output times in the range 
G.17- 13.5 Gyr 



Following iLacev &: Cold il993l) we define the formation 
time of a halo as the time when it has assembled half of 
its final mass. The distribution of formation redshifts for 
a halo with final mass Mi at redshift zi can therefore be 
approximated by 



(7) 



p{z{)dz{ — p{Lo)dLu — 2lo erfc [ —jzz ) do). 



where lo'^ = {Sd - (5ci)^/(Sf - Si), (5cf = &c{zi), 5^ = Sc{zi) 
and Si — cr^(Mi/2). Here Sc{z) is the overdensity required 
for spherical collapse at redshift z extrapolated using linear 
theory to the present time and a^{M) is the variance of the 
initial density fiuctuation field. Using this description we can 
estimate (approximately) the age distribution of haloes as a 
function of Wvir at each output redshift. At high redshifts the 
halo age distribution is of course narrow and broadens sig- 
nificantly at lower redshifts where a large number of output 
ages are needed to model the halo population. 

At each output redshift we calculated the number den- 
sity of haloes as a function of virial velocity using the 
ISheth fc TormenI (|T999fl mass function 



nsT{M)dM = - 



(I) 



1/2 



1 + 



1/2 Pb 
M 



Sc da 
^ -TT7 exp 
cr2 dM 



a5l 

'2C72 



dM 



(8) 



where A = 0.322, p = 0.3, a = G.707, pb is the background 
density, a is the rms density fiuctuations on scale R corre- 
sponding to mass M (M — 47rp6_R"^/3) and 5^ ~ 1.68. 

The distribution of halo spin parameters is assumed to 
follow a lognormal distribution and is assumed to be inde- 
pendent of both time and Uvir 



p{X)d\ ■ 



■ exp 



log' (A/A) 



2al 



dX 
T 



(9) 



with A = G.045 and ax = 0.56 ijVitvitska et alJl20G2h . 

All haloes are assumed to evolve in isolation. We have 
not included an elaborate merger tree to describe the merger 
history of an individual dark matter halo and its baryonic 
content. Our model is thus not designed to give a descrip- 
tion of the history of an individual galaxy. Instead our model 
has been constructed to provide an approximate 'snapshot' 
of the galaxy distribution at each output redshift character- 
ized by the parameters (tf,Wvir,A). By using the simplified 
model of equations ^ and JSJ we are, in effect, assuming 
that the baryonic parts of galaxies follow a similar merger 
history to that of the dark matter within the virial radius. 
This is, of course, a simplification and one would expect the 
high density baryonic components to survive as distinct sys- 
tems for a longer time than their parent haloes. Our main 
results on the metallicity distributions of DLAs (Section 3) 
are insensitive to the merger history. However, some prop- 
erties, in particular the numbers of DLAs at high redshift, 
are extremely sensitive to the merger history (see Section 4). 
Our model can give only a very approximate description of 
properties that are sensitive to merger history and further 
progress will probably require high resolution gas dynamical 
numerical simulations. 

At each output redshift we calculate properties of the 
galaxy distribution along sight-lines assuming galaxies are 
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Figure 2. The ejection fractions defined as /q = Afcj/Mtot, where Mtot = Mgj + Mgas + A/star is the total baryonic mass. The results 
are plotted at z = (left panel) and z = 3 (right panel) as a function of the spin parameter A and the virial velocity of the halo v^i,. for 
the 119 models. The ejection fraction shows a strong increase with decreasing A and Dvir- 
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Figure 3. The metallicity of the ejected gas relative to the galaxy gas metallicity (Z^j/Zg) at z = (left panel) and z = 3 (right panel) 
as a function of the spin parameter A and the virial velocity of the halo Dvir for the 119 models. The ejected gas is predominantly enriched 
relative to the gas in the parent galaxy. 



uniformly distributed and have random inclination angles. 
We assumed a hydrogen fraction of X = 0.7 by mass. If the 
probing sight-line encounters a hydrogen column density of 
A''hi ^ 2 X 10^°cm~^ « 1.6 M0pc~^ we select the sightline as 
a DLA and measure its physical properties such as column 
density, metallicity and star formation rate at the intercept 
point. 



2.3 Model properties 

In our model, galaxies form in an intense starburst followed 
by a rapid buildup of metallicity, as can be seen in Fig. 
which shows the formation of a disc in halo with a relative 



high A and a virial velocity Vvir = 100 kms~^ (see also EOO 
for a more extensive discussion) . The initial burst is followed 
by quiescent star formation that builds up the galactic disc 
from inside-out. The model develops an extended gaseous 
disc with surface density well above the DLA definition of 
A'hi 2 X lO^'^cm"^. The sharp outer edge of the gas disc 
seen in Fig. is a consequence of the infall model in which 
the final time of the model sets the maximum cooling ra- 
dius within the halo. The oscillatory behaviour of the radial 
gas profile near the centre is a consequence of the galactic 
fountain and similar behaviour can be also seen in the star 
formation profile. 

In addition to gaseous infall and star formation, the 
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evolution of galaxies according to this model is strongly 
influenced by outflows produced by stellar feedback. The 
ejected mass fraction, /oj, at time t, is defined as the ra- 
tio of the ejected baryonic mass Moj divided by the total 
baryonic mass Mtot = Moj + Mgas + Mstar, where Mgas and 
Mstar are the mass in the cold gaseous disc and the stel- 
lar disc at time t. The parameter /ej is a strong function of 
both Vvii and A as can be seen in Fig. |5| According to our 
model, galaxies with small v^ir and A can have 75% or more 
of their gas mass ejected by 2; = 0. This is expected to have a 
large effect on various properties (particularly metallicities) 
of DLAs as the low mass haloes dominate the cross-section 
especially at high redshift (see Fig.|5J. 

In addition to their large ejection fraction, galaxies with 
small Uvir and A predominantly eject gas that is enriched 
relative to the host galaxy. This is especially true at high 
redshifts as can be seen in the right-hand plot of Fig. |21 By 
z — the metallicity enhancement of the ejected material 
is lower because most of the ejection occurs at high redshift 
and the galactic disc, by z = 0, has had time to evolve to 
higher metallicity. 

The effectiveness of supernova feedback in ejecting gas 
from spiral and irregul ar galaxies has also been studied ob- 
servationally bv iGarne tt (2002). They found an increasing 
correlation for both the metallicity Z and the effective yield 
yeff as a function of rotational velocity Vyot- However, at 
Vrot ~ 150 km/s both correlations turn over and both the Z 
and UeS becom e approximately constant for v^ot ~ 150 km/s. 
iGarnetj (|2Qq3) conclude that this observed trend suggests 
that galaxies with Wrot ~ 100 — 150 km/s may lose a large 
fraction of their metal-enriched gas through supernova feed- 
back, whereas more massive galaxies tend to retain their 
metals. 

To test this observational picture we calculated the av- 
erage stellar metallicities of our model galaxies as a func- 
tion of Uvir and A at different redshifts. Our models pro- 
duce an increasing correlation of the stellar metallicity with 
Vvii until a critical Vvii is reached, after which the metal- 
licity becomes approximately constant. At z — this crit- 
ical velocity is Uvir ~ 70 — 85 km/s, which corresponds to 
Vrot ^ 170 — 200 km/s for a typical model. Thus the chem- 
ical evolution of the stellar component in our model seems 
to be in very good agreement with observations. 



3 METALLICITY EVOLUTION OF DLAS 

3.1 Normal feedback (NFB) model 

The metallicity evolution of the DLAs as a function of red- 
shift is shown in Fig. |1] for the normal feedback model 
with (j>K_ = 0.1. The square (green) points show the mean 
DLA metallicity and the (blue) shaded contours illustrate 
the metallicity distributions at each output redshift. The 
diamond (yellow) points show the predicted A^Hi-weighted 
mean metallicities. The metallicity evolution binned in unit 
redshift bins for the observational data set and our feedback 
model predictions is summarized in Table I^TI 

Figure |1] shows observational measurements of [Zn/H] 



for DLA systems^. The observational sam ple consists of 
87 ab undance measurements compiled by iKulkarni et alJ 
I2OO5I) (see ref e rences therein) and 16 measurements by 
lAkerman et al.l ^2005^ totalling 103 abundance measure- 
ments. These data consist of both observations with error 
bars and upper limits indicated by triangles. We used the 
[Zn/H] abundance as an indicator of total metallicity to 
compare with our models. The assumption that zinc is un- 
depleted on grains and traces the iron a bundance over three 
orders of magnitude in [Fe/H] (see e.g. lSnederi et alJll99J) 
is stro ngly supported by the recent work of l^^senet^D 
(l2004f) . Zinc should therefore provide an accurate overall in- 
dicator of the total metallicity in the gas. In addition, by 
using [Zn/H] abundances we can compare our results to a 
large observational sample. 

The most important result from Fig.|3]is that our model 
predicts a typical DLA metallicity of [M/H] « —1.0 with 
very little evolution over the entire redshift range z — 5 
to z — 0. In more detail, the model predicts a mean A^^hi- 
weighted metallicity of [M/H] — —0.86 at z = 5 and shows 
a slight decrease from z — 5 to z ^ 2. This decrease can be 
seen in the model shown in Fig.Qand is accentuated by the 
stronger outflows of enriched gas from galaxies with small 
Hvir ~ 70 km/s which dominate the DLA cross-section at 
high redshift. From z — 2 to z = the model shows a small 
but steady increase in the average metallicity. This is caused 
by two effects: a) galaxies at lower redshifts have had time to 
build up higher metallicities via quiescent star formation; b) 
the DLA cross-section becomes dominated by galaxies with 
larger virial velocities («vir = 70 — 100 km/s) for which the 
effects of metal enriched outflows are not as pronounced as 
for galaxies with lower Uvir. The combined effect of a shift in 
the typical v^ir of the DLA population with redshift and the 
dependence of the metallicity of the outflows as a function 
of Vvir and time results in almost no evolution in the DLA 
metallicity. It is worth emphasising that this follows with no 
adjustment of the feedback parameters compared to those 
used by EOO. The parameters of the model have not been 
tuned specifically to match the observations plotted in Fig. 
13] instead we adopted a model that roughly matches the 
gross properties of the ISM of the Milky Way and reproduces 
its net star formation rate. 

Performing a least-squares linear fit to the mean A^hi- 
weighted model predictions over the entire redshift range 
we derive a slope of m = —0.04 dex/Az with a zero point 
of [M/H]o = —0.76 dex confirming the visual impression 
of no evolution within the error bars. The fit to the un- 
weighted mean metallicities yields m = —0.02 dex/Az with 
[M/H]o = —1.26 dex. The unweighted metallicity tracks the 
Ani-weighted metallicity accurately with an offset of 0.5 
dex. 

However the model predicts some evolution in the 
metallicity distributions of the DLAs with a tail extending to 
low metallicities at redshifts z ^ 1.5. In fact our model pre- 
dicts that about 5%-10% of the DLAs should have metallici- 
ties of [M/H] < —3 in the redshift range z — 2—4. In contrast 
there is an apparent observed 'floor' in D LA metallicities. 
In particular, in the sample analysed by IProchaska et alJ 



The observational data points were kindly provided by Chris 
Akerman 
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Figure 4. The redshift evolution of DLA metallicity for our normal feedback model (</)k = 0.1) overplotted with observations from 
I Akcrman ot al. 2005; Kulkarni et al. 2005) and references therein. The (blue) shaded contours give the range of metallicity for DLAs 
at each output redshift. The square (green) points show the mean DLA metallicity and diamond (yellow) points show the Afpji-weighted 
mean metallicity. 
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Figure 5. The redshift evolution of DLA metallicity for our weak (left-hand panel) feedback model (</<« = 0.01) and strong (right-hand 
panel) feedback model {4>k = 1-0) overplotted with the same observational data and symbols as Fig. El 



l|2003h no single DLA system was observed with a metallic- 
ity [M/H] ^ — 3 . These authors use a-elements and Fe mea- 
surements, in addition to Zn measurements, and so are able 
to detect lower metallicities than would be possible with Zn 
alone. At face value this suggests a discrepancy, however our 
metallicity model assumes that the primordial infalling gas 
has zero metallicity. By allowing for a mild pre-enrichment 
of the intergalactic medium to the level of [M/H] = —4.5 



the 99%-contour would move to [M/H] ~ -4 with the 95%- 
contour also at higher metallicity resulting in better agree- 
ment with the apparent observed metallicity floor. We will 
return to this issue in Section 3.3, where we discuss possible 
effects of selection biases and dust obscuration. 

Comparing to the observations, we see that the majority 
of the data points lie within the 75%-contour with some out- 
liers at high metallicities. Especially at around z — 2 a few 
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Data set < z < 0.5 0.5 < z < 1.5 1.5 < z < 2.5 2.5 < z < 3.5 

AUobs -0.87 ±0.25 -1.02 ±0.27 -1.00 ±0.15 -1.20 ±0.18 

Errbars -1.06 ±0.28 -1.01 ±0.15 -1.17 ±0.28 

NFB -0.70 -0.88 -0.90 -0.90 

WFB -0.66 -0.70 -0.66 -0.60 

SFB -0.77 -0.99 -1.13 -1.19 



Table 1. The metallicities as a function of redshift for the observations and for our three feedback model predictions. 'All obs' includes 
all observations, with upper limits as if they were detections, whereas 'Errbars' excludes the upper limit data points. The observations 
have been included with 95% uncertainties given by a bootstrap analysis. For the lowest redshift bin we only have one observation with 
errorbars. 
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Figure 6. The HI column density versus gas metallicity at redshifts z=0-5 for the NFB model (left-hand panel) marked with N and for 
the NFBc coeval evolution model (right-hand panel) mar ked with C. The model (open circles) is over plotted with observations with error 
bars (red squares) and upperlimits (red triangles) from jAkerma n et_alj|200^jKulkarn i et alj|2005l) and references therein. In addition 
we include the recent gamma-ray burst DLA (red star) observed~bv lFvnb^e^!Ln20o3) . The observational data has been binned in unit 
redshift bins at each redshift. 



DLAs have been observed with 'anomalously' high metallici- 
ties of almost sol ar value. In fact, the h ighest metallicity sys- 
tem observed bv lLedoux et al.l (|2002|) has [M/H] = -0.11. 
These authors argue that this DLA is probably connected 
to a star formation region because they find evidence of 
dust and a relatively high molecular hydrogen fraction. Our 
model does not include spatial inhomogeneities within galax- 
ies which could affect the predicted metallicity distributions. 
Furthermore our model does not include any variance in the 
strength of stellar feedback, which as we will show in the 
next subsection could also influence the predicted distribu- 
tions. 

In summary the results of this Section show that our 
model reproduces the lack of evolution in the metallicity of 
the DLA systems seen in the observations over the wide red- 
shift range z « 4 to z = 0. In addition, the mean metallicity 
and dispersion predicted by the model also provides a good 
match to the observations. There may be some discrepan- 



cies between the predicted metallicity distributions and the 
observations at both extremely high and low metallicities. 
However as we will argue later in this Section the tails of 
these distributions are sensitive to model parameters and 
selection biases. 



3.2 Dependence on the strength of stellar 
feedback 

We also calculated the metallicity evolution adopting strong 
{<f>K ~ 1.0) and weak (0„ = 0.01) feedback prescriptions. 
The predicted metallicity distributions as a function of red- 
shift are plotted in Fig. |K| As expected, the SFB model 
produces strong outflows that eject a high fraction of the 
metal enriched gas from the galaxies, especially in their 
early phases of evolution. The overall effect is to produce 
a net metallicity that is 0.1-0.3 dex lower than for the NFB 
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Figure 7. The redshift evolution of DLA metallicity for the bias-corrected normal feedback model (cp^ = 0.1) overplotted with the same 
observational data and symbols as Fig. El 



Data set 


m 


[M/H]o 


NFB 


-0.04 


-0.76 


NFB unweighted 


-0.02 


-1.26 


WFB 


-1-0.03 


-0.69 


SFB 


-0.07 


-0.88 


NFB bias-corrected 


-0.03 


-0.79 


NFB dust-corrected 


-0.03 


-0.84 



Table 2. The best fit linear fits to the metallicity evolution 
with redshift for the various data sets, m is the slope in units 
of dex/Az and [M/H]o the zero point in units of dex. All metal- 
licities weighted with Nm, unless otherwise stated. 



model, with the difference being larger at high redshifts. 
Fitting the linear evolution of the A'ni-weighted metallicity 
we get m = -0.07 dex/Az with [M/H]o = -0.88 dex. This 
model produces slightly stronger evolution in the metallicity 
because stronger outflows lower the metallicities at high red- 
shift in comparison to the NFB model. However, the main 
effect of the strong feedback is to lower the overall mean 
metallicity, making it more difficult to explain the observed 
high metallicity DLAs. 

The WFB model produces mean metallicities that are 
0.1-0.3 dex higher than for the NFB model. The main dif- 
ferences with the NFB model are at high redshift where 
the reduced effects of outffows allow the buildup of high 
mean weighted metallicities of [M/H] ~ —0.5. This high 
early enrichment results in a small positive evolution in 



metallicity with increasing redshift m — -1-0.03 dex/Az with 
[M/H]o = -0.69 dex. 

Both the strong and weak feedback models produce es- 
sentially no evolution in the metallicity. The main effect of 
altering the feedback efficiency is to raise or lower the over- 
all metallicity by a relatively small amount. Interestingly 
our default NFB model comes closest to providing a match 
to the observations. Evidently the assumption of a uniform 
feedback efficiency is a gross oversimplification. More real- 
istically, we would expect a range of feedback efficiencies 
caused by both differences in the properties of the interstel- 
lar medium in galaxies (cold cloud temperatures, sizes, heat- 
ing mechanisms etc.) and from transient departures from our 
self-regulated model of star formation {e.g. stronger stellar 
feedback in interacting systems and weaker feedback in more 
isolated systems). Variations in feedback efficiencies might 
well affect the overall metallicity distributions of DLA sys- 
tems particularly in the tails at high and low metallicities. 
However, although it would be possible to construct an ad- 
hoc model with variable feedback efficiencies, developing a 
physically motivated model is beyond the scope of this pa- 
per. 



3.3 Selection biases and dust 

The model predictions of the previous sub-section were com- 
puted assuming that the eff ects of dust and oth er observa- 
tional biases are negligible. iBoisse et alJ lll99^ report an 
anti-correlation between the observed Zn abundance and 
A'hi, independent of redshift. This result has been confirmed 
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by IPrantzos fc Boissierl i200Ctl who note that all observed 
DLAs are found between 18.8 < [Zn/H] + log(NHi) < 21 
in the [Zn/H] — Nhi plane. The lower limit in this relation 
may arise because low column density and low-metallicity 
systems are not detectable through their absorption lines 
in current surveys. There are several possible explanations 
for the lack of DLAs with high column density and high- 
metallicity: a) The light of the background quasars may 
be severely extinguished by dust; b) the cross-section for 
high A'^Hi absorption may decrease with increasing metallic- 
ity as gas is consumed by star formatio n; c) system s with 
log(A^Hi) > 21 may be intrinsically rare JSchavdl200ih . 

The left-hand panel of Fig.|S|shows our default (j)^. = 0.1 
model predicti ons in the [Zn/H] — Nhi pl ane. The straight 
lines show the iPrantzos fc Boissierl ll2000ll bias thresholds. 
The fraction of our models that lie above the upper cutoff 
(and therefore may be biased by dust obscuration) is always 
low except for the higher redshifts where its contribution 
is of order 7% — 8% ai z — 4 — 5. However we find that a 
significant fraction of our models lie below the lower thresh- 
old (typically 10-20%). It is therefore clear that removing 
points below the lower line will have a significant effect on 
metallicity distributions, whereas removing points above the 
upper line will have little effect. We have also overplotted 
the observational data on this diagram. At redshifts z < 1 
the observations match our models reasonably well but at 
z > 2 the models fail to match the observational points 
with logA^Hi < 20.5 and metallicities [Zn/H] > —1. In the 
bin centered on z = 4 there is almost no overlap between the 
observations and our NFB model, though at high redshift all 
but one of the quasar selected metallicity observations are 
upper limits. We have also plotted the recent metallicity de- 
termination for a DLA system ai z = 4.048 detected against 
the optical after-gl ow of the gamma-ray burst GRB060206 
jFvnbo et alJl2006D . This lies much closer to the model pre- 
dictions. 

Recently there has been consider ably interest in sub- 
damped Lyman-Q systern s (see e.g. iLedoux et alJ 120031 : 
iDessauees-Zavadskv et all 12003') with column densities of 
lO^^cm"^ < Nhi ^ 2 X 10^°cm"^, i .e. below the classi- 
cal definition of DLAs. The analysis bv lPeroux et alJ ll2005ll 
found that the sub-DLAs contribute about 20% to the total 
neutr a,l gas mass de i isity at redshifts 2 < z < 5. Further- 
more, iPeroux et alJ i200d) recently discovered a sub-DLA 
with super-solar metallicity at z = 0.716 indicating that sub- 
DLAs might be more metal-rich than classical DLAs. We do 
not consider sub-DLAs in the present analysis. However, we 
note that a large number of sub-DLAs with low-column den- 
sities and high metallicities cannot be explained within the 
framework of our model as discussed below. 

Our model cannot produce large numbers of low col- 
umn density systems at high redshift with high metallicity. 
As can be seen from Fig. the model tends to produce 
flat A'hi distributions of high column density with a sharp 
cutoff. Clearly, the sharpness of the cut-offs are an ideal- 
ization and a consequence of applying a simplistic model 
of self-regulated star formation based on the stability cri- 
terion of equation Q. In reality, there are many physical 
process that will tend to 'smear' the cut-offs. These include 
tidal disruption and realistic deviations from axial symme- 
try both during collapse, and from disc instabilities. (In this 
context, it is worth noting that the highest column densi- 



ties within the Magellanic Stre am are only just belo w the 
DLA column density threshold jPutman et al. ll2003ll . sug- 
gesting that tidal debris may contribute to the DLA cross 
section). Additional effects, not incorporated in our model, 
that might increase the spread in column densities at high 
metallicities include bursts of star formation (triggered by, 
say, an interaction) that can consume much of the gaseous 
disc, or ram pressure stri pping of disc gas. We not e that the 
numerical simulations of 'Naga mine et aP |2004b) produce 
a wide spread of column densities extending down to the 
DLA threshold. However, they find a strong correlation be- 
tween metallicity and column density and (as in our models) 
cannot account for the low metallicity, high column density, 
systems observed at high redshift. Understanding the col- 
umn density distributions plotted in Fig.|S]in greater detail 
poses an interesting challenge for both theory and observa- 
tions. 

Our simplified model for the age distribution of disc 
galaxies can affect some aspects of our model at high red- 
shift. This is discussed in greater detail in Section 4, but to 
illustrate the effects on the column density distributions we 
plot in Fig. |S| the evolution of a model with coeval evolu- 
tion and a formation redshift set to infinity (model NFBc) 
in add ition to our standard models using the lLacev fc Cola 
(Il993h (hereafter LC93) age distributions. In the coeval 
model, discs have more time to evolve, especially at high 
redshifts as discussed in Section 2] The additional evolution 
time shifts the A'^hi distributions towards lower column den- 
sities, (and the gas metallicities to slightly lower values) but 
not by enough to match the observations at z ^ 3. 

In Fig.|7|we plot the metallicity evolution with redshift 
for the DLAs in our models that lie within the diagonal 
lines plotted in Fig. El i.e . the relations that according to 
IPrantzos fc BoissieJ (l2000h roughly delineate regions in the 
diagram that might be affected by observational selection 
biases. Including the lower limit, 18.8 < [Zn/H] + log(NHi), 
makes quite a large difference because it eliminates most 
of the systems with metallicities less than [M/H] — —2. 
After this bias-correction almost all of the observations at 
low metallicity lie within the 99% contour predicted by 
our model. It may therefore be necessary to understand a 
bias of this type in more detail to interpret the Zn abun- 
dances at very low metallicities. This type of bias will be 
less important in direct observations of Fe, which is more 
abundant th an Zn. Indeed a s we have mentioned previ- 
ously tProchaska et all l2003il find abundances as low as 
[M/H] ~ — 3 using Fe and a elements. 

The bias-correction at high column densities removes 
some DLA systems in our models with high metallicity, 
which makes it even more difficult to explain the observed 
high metallicity systems at z ~ 2. Fitting a linear evolution 
model to the 'bias-corrected' model predictions of Fig.jTjre- 
sults in m = -0.03dex/Az with [M/H]o = -0.79 dex, i.e. 
virtually identical to fitting the uncorrected results shown 
in Fig. El 

If we only impose the upper cut, [Zn/H]-|-log(NHi) > 21, 
the fit to the metallicity evolution gives m = —0.03 dex/Az 
with [M/H]o = —0.84 dex. In this case our model predicts 
a slightly lower metallicity of ~ 0.1 dex compared to the 
uncorrected results of Fig. 21 According to our model only 
a small number of systems are expected to lie above the 
line defined by the line [Zn/H] -|- log (Nhi) > 21, and re- 
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moving the ones that do lie above this line has very little 
effect on the overall metallicity distribution. Dust obscu- 
ration in high metallicity, high column density systems, is 
therefore unlikely to lead to significant biases in the DLA 
metallicity distributions according to our model, since such 
systems are ra r e. Th is is consistent with the results of 
lAkerman et al.l (|200^, who find only a marginally higher 
metallicity of 0.2 dex for their unobscured (radio-selected 
quasar) CO RALS DLA sampl e compared to a control sam- 
ple from iKulkarni et alJ (l2005l) . 



4 GLOBAL PROPERTIES OF DLAS 

In this Section we compute various global properties of the 
DLA population. Many of the properties discussed in this 
Section are sensitive to the merger history of the absorbers, 
and hence to the age of the model galaxy that determines 
its cross-section and mass. As explained in Section 2.2, we 
have employed a simplified description based on the LC93 
model for the age distributions of dark matter haloes as a 
function of their virial velocity. The high density baryonic 
systems within these haloes will be longer lived than the 
haloes themselves and will follow a more complex merger 
history. The metallicity distributions described in the previ- 
ous Section are insensitive to models of the merger history 
because the metallicity of the gas is not a strong function of 
age. Likewise, the gas metallicity gradients in our models are 
shallow and so the metallicity distributions are insensitive to 
any process, such as tidal stripping, that might truncate the 
gaseous discs. The metallicity distributions discussed in the 
previous Section are therefore robust and are largely a con- 
sequence of the model for feedback and self-regulated star 
formation. 

The cross-sections, on the other hand, are extremely 
sensitive to the age distribution and hence to the merger 
history. Although it is possible, in principle, to construct a 
more elaborate a nalytic merger model for the b aryonic com- 
ponents (see e.g. IPenarrubia and BensorJl200.4l the physics 
involved is complex and difficult to model accurately. To 
illustrate the sensitivity of various results to the merger his- 
tory, in addition to our standard models using the LC93 age 
distributions, we have therefore calculated the evolution of a 
model with coeval evolution and a formation redshift set to 
infinity (model NFBc) . The NFB and NFBc models provide 
two extreme representations of the merger history of early 
disc systems that should bracket the predictions of a more 
realistic model. The differences between the two represen- 
tations is particularly important at high redshift. At z — 5 
the age of the coeval population is 1.2 Gyrs compared to the 
LC93 model which predicts ages for the models of only 0.3 
Gyrs. This difference of a factor of four leads to large dif- 
ferences in the calculation of cross-sections, redshift number 
densities and Qm. 

To summarize, our model provides robust predictions 
for the metallicity evolution of DLAs. Properties that are 
sensitive to the disc age distributions and cross-sections, 
particularly the rates-of-incidence (see equation llOfl and 
(to a lesser extent) fini are less reliably predicted by our 
model, especially at high redshifts. Many of the results pre- 
sented in this Section are therefore meant to be indica- 
tive at the 'order-of-magnitude' level, rather than detailed 



fits to the observations. To make accurate predictions of 
such properties, it is likely that high resolution numerical 
simulations that properly model the formation, evolution 
and merger history of dwarf galaxies will be required ( cf. 
[Naaamine et al...2004al . 



4.1 Abundance of DLAs as a function of column 
density and redshift 

Follo wing iNagamine et alJ (|22^3) (see also I Wolfe et'ai] 
I2OO5I and references therein) we can compute the total DLA 
rate-of-incidence from 



diVoLA, 

dz 



(2) = ^/ ndn.(M',z)c7DLA(M',z)dM', (10) 



where ndm(W,z) is the ISheth fc Tornienl lll999l) (comov- 
ing) dark matter halo mass function, (Tdla(M',2;) is the 
comoving DLA cross-section and dr/dz — c/H{z), with 
H{z) = Ho^/nmil + z)s + nA. The results for our models 
are plotted i n the left panel of Fig. El together with obser- 
vations from ( Prochaska fc Herbert-For3l2004l : IZwaan et alJ 
12005b.:, Rao et ah. 2006') . 

We discuss first the results for our standard models 
based on the LC93 age distributions. The strength of feed- 
back has almost no effect on the rate-of-incidence, but the 
overall shape is very different from the observations. At red- 
shifts below two the models match the observations reason- 
ably well (within a factor of two) but at higher redshifts, the 
models significantly underproduce dNoLA/dz by a factor of 
five. The main reason for the underproduction of dNuLA/dz 
at high redshifts is that the prescription of LC93 leads 
to galaxies of very young ages (and therefore small cross- 
sections) at high redshifts. Assuming coeval evolution as de- 
scribed above (model NFBc) produces dNuLA/dz which is 
higher by a factor of five at high redshifts, roughly match- 
ing the observational data at z ^ 4. However, the coevo- 
lution model significantly overpredicts dNuLA/dz at lower 
redshifts. 

As we have argued above, the true age range of the 
galaxy distribution is likely to lie somewhere between the 
LC93 and coeval distributions and so Fig.|Hlsuggests that it 
may be possible to account for the observed rate-of-incidence 
with a model of the type described here incorporating a 
m ore realistic merger his tory. The numerical simulations 
of INagamine et all ^2004a^ actually reproduce the observed 
dNY)-LA/ dz quite well and in fact their comoving tabulated 
cross-sections as a function of Wvir lie in between the cross- 
sections of our NFB and NFBc models. It would be interest- 
ing to make a more detailed comparison of our models with 
age distributions and cross-sections from numerical simula- 
tions to determine the key physics required to reproduce the 
observed rate-of-incidence. 

The right-hand panel of Fig. |H] shows the differen- 
tial distribution of f{NYa,X{z)) = d^N/dNmdX pre- 
dicted by our models. This quantity gives the number 
of DLAs per unit column density A ^ hi an d per unit 
absorption distance ('Bahcall fc Peeble^ ^6^. The data 
points plotted on the righ t -hand panel of Fig. |H| are from 
IStorrie-Lombardi fc WolJM ll2000lh adjusted for our ACDM 
cosmology. The model results are calculated by averaging 
f{Nm) over the redshifts ^ z ^ 5 and all three models 
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Figure 8. The DLA rat e-of-incidence dNY)i^\{z) / dz (left- hand panel) for the thre e feedback model s and the coeval evolution model 
is plotted with data from |Prochasla_& Herbe^ l2004^ (solid circles at 2 > 2'). iRao et alj l2006l) (solid squares at 2 ^ 1) and the 

open square at 2 = froni Tzwaai^^aLrT2005bF ~Tri the right p anel the DLA column density dist ribution /(A'^hi) is shown for the four 
models overplotted with observational data (solid symbols) from lStorrie-Lombardi fc Wolf j l200Cl) . The model results are averaged over 
^ 2 ^ 5. 



fit the observational data reasonably well at high column 
densities but underproduce systems at the low column den- 
sity end by a factor of two. The three models differ only 
at the highest column densities, where the weak feedback 
model produces slightly more high column density systems 
compared to the normal and strong feedback models. On the 
other hand, the coeval evolution model fits the observational 
data nicely with the exception of the highest column density 
point. Again, as with the rate-of-incidence analysis, the two 
extreme merger models roughly bracket the observational 
data. 



4.2 DLA halo properties 

In the two next sub-sections we discuss the properties of 
typical DLA haloes and of the galaxies that lie at their cen- 
tres. Unless otherwise stated, we show results for NFB model 
using the LC93 halo age distributions. 

The evolution of the average virial velocity of the DLA 
haloes, together with the 50% (quartile, black lines) and 
95% (red lines) ranges, are shown in Fig. We find 
a gradual but mild increase of the average velocity from 
■Uvir ~ 60 kms^^ to Uvir ~ 75 kms~^ with a large scatter. 
The contribution of systems with Hvir ^ 100 kms~^ to the 
DLA cross-section increases from 14% at z = 5 to 22% by 
2 = 0. The contribution of even higher virial velocities of 
VviT ^ 200 kms~^ (i.e. haloes hosting large Milky- Way type 
galaxies) remains low at less than 2% at all redshifts. Our 
re sults agree we l l with those from the semi-analytic models 
of lOkoshi et alj (|20o3), especially at high redshift. At lower 
redshift, the typical virial velocities are somewhat lower than 
their s. However our result s for z ~ are in better agreement 
with IZwaan et alJ ^2005b^ who infer a mean virial velocity 
for the local DLA population of Vviv ~ 70 kms~^. 

In Fig. 1^ we show the corresponding evolution in the 



median mass (within the virial radius) of the dark matter 
haloes containing DLAs. We find that the median mass in- 
creases from Mdm ~ 5 x l(f Mq at 2 = 4 — 5 to a few times 
1O^°M0 at « = 0.5 - 3.5 with an upturn to Mdm ~ lO^M© 
at 2 = 0. Our halo masses are somewhat lower than those 
inferred by iBouche et alJ (jioO^) in an observational study 
of Mgll absorbers in the redshift range 0.4 ^ 2 ^ 0.8. 
Their results imply that the hosts of Mgll absorbers (typ- 
ically 40% of Mgll absorbers are DLAs) have halo masses 
of ~ 2 — 8 X 10^^ M0, perhaps suggesting that Mgll selec- 
tion biases samples towards higher virial velocities than the 
average. The abrupt increase at 2 ^ 0.25 is caused by the 
exhaustion of gas in low mass haloes which placing them 
below the DLA selection criterion of A^^hi ^ 2 x lO^^cm"^. 

We conclude that the DLA cross-section is dominated 
at all redshifts by systems residing in small mass (low virial 
velocity) haloes with a mild evolution towards systems with 
higher mass at lo w redshift. This is in agreement with re- 
cent semi analvtic | CTcoshi et alJl2004r) and numerical inves- 
tigations llCen et ahlbOOal) but in strong disagreement with 
the classical picture in which the majority of DLAs at all 
redsh ifts are massive disc galaxies like our own iWolfe et alJ 
Il986l) . 

The spin parameter A of the DLA haloes as a func- 
tion of redshift is shown in Fig. |^. As DLAs are selected 
by cross-section, one would expect DLAs to be biased to- 
wards bigger discs with larger A in comparison to the gen- 
eral galaxy population. This is seen in our models. The DLA 
cross-section is dominated by systems in haloes with large 
A, especially at high redshift where high A systems rapidly 
build up large gaseous discs. Systems with small A also ex- 
perience more extensive mass loss caused by the feedback 
from supernovae (see Fig. At lower redshifts the mean 
cross-section weighted A decreases. This is explained by the 
slower buildup of gaseous discs in lower A systems. In sum- 



© 2005 RAS, MNRAS OOO.imTOl 



A model for the metallicity evolution of damped Lyman-a systems 13 



200 



\ 150 - 
B 

A 

100 - 

V 



50 - 



0.15 - 




1 2 3 4 5 
z 



1 2 3 4 5 
z 



O 





Figure 9. The evolution of the mean v^ir and A with 50% (quartile, black lines) and 95% (red lines) ranges (top-panel). Evolution of 
the median DLA halo mass and comoving halo density with 50% (quartile, black lines) and 95% (red lines) ranges (bottom-panel). For 
points where only one low/high error bar is shown the 50% and 95% ranges are the same, which is a result of the limited number of bins 
in our model (17 in Hvin -'^DM a^nd ™dla; 7 in A, see Section 1^.21 for details). All plots arc for the normal feedback model (0^ = 0.1). 



mary, we conclude that the DLA cross-section is dominated 
by systems with larger spin parameters than the mean, and 
that the bias is accentuated at high redshifts. 



The comoving number density of the DLA parent haloes 
is shown in Fig. |^. Mirroring the evolution of the typical 
DLA halo mass with tidm ~ 0.1 Mpc"'' at z = 4 — 5 de- 
creasing to riDM ~ 0.01 Mpc~^ by z=0. This shows that 
the DLA population acquires an increasing contribution at 
lower redshifts from more massive haloes with lower comov- 
ing number densities. We can see a jump in both the masses 
and number densities between z = 5 and 2 = 4, this is again 
caused by the slower buildup of the smallest mass systems, 
which only start to fully contribute to the DLA cross-section 
at z ^ 4. 



4.3 DLA galaxy properties 

Fig. llOl shows the evolution of various properties of the DLA 
galaxy population. The DLA median hydrogen gas mass 
evolves by a factor of ten from Mm ~ IO^Mq at 2 = 5 to 
Mhi = 2 X 10^ Mq at z = 0. Our model prediction at 2 = is 
in excellent agreement with lZwaan et all i2005lJ) . The evo- 
lution of AIhi is associated with an increase in the median 
stellar mass from Af* ~ 6 x W"^ Mq at z = 5 to M* ~ 10^ 
by z = 0. According to our model, the DLA cross-section is 
dominated by 'dwarf (low-mass) galaxies with large HI discs 
at all redshifts and especially at higher redshifts. However, 
the width of the distribution (indicated in the figure) in- 
creases at lower redshifts because higher mass systems make 
a growing contribution to the DLA cross-section. The con- 
tribution of Milky Way-type systems, with total disc masses 
of Af ^ 5 X lO^^AfQ, to the cross-section is completely neg- 
ligible at high redshifts, growing to only 0.5% at z = 2. At 
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lower redshifts this fraction grows to 2% at z = 1 with an 
increase to 6% by the present epoch. These result s agree 
well with the semi-analytic models of lOkoshi et alJ i2004l . 
who found a ~ 5% — 10% contribution of Milky Way type 
systems at low redshift with negligible contribution to the 
cross-section at z ^ 3. According to our model, the major- 
ity of the DLA population must consist of 'dwarf galaxies 
with typical masses of ^ O.IM*. This result is also consis- 
tent with the lack of direct detections (see e.g. IWolfe et alJ 

Figure IlUt shows the evolution of the median impact 
parameter (in physical units) as a function of redshift. We 
find a very strong, almost linear evolution in this parameter, 
from very small values of ~ 1 kpc at z = 5 to values of 
~ 7 kpc at the present epoch. The local value again agrees 
well with the observational results of .Zwaa n et al. (2005bi) 
who find a median value for the local impact parameter of 
7.8 kpc. The small impact parameters in combination with 
the relatively low number densities of the DLA haloes at 



high redshifts results in the low rate-of-incidence at high z 
discussed in Section 4.1 (Fig. (HJ. 



4.4 Star formation of DLAs 

The typical star formation rates of DLAs according to our 
model are shown in Fig. IIUH . We find that the DLAs have a 
very broad distribution of star formation rates, ranging from 
~ 5 X 10"^ up to ~ 50 Meyr"^ at all redshifts. The evolu- 
tion of the median star formation rate is virtually constant 
with a median star formation rate of Af* ~ 0.1 — 0.2A'/Qyr~^ 
over the entire redshift range z = 0-5. The mean star forma- 
tion rate is considerably higher at M* ~ 1 — 2Af0yr~^. The 
distribution of star formation rates has a long tail extending 
to high values for sightlines probing the central regions of 
massive galaxies (see the 95% ranges in Fig. llOti ). 

However, we find strong evolution in the star forma- 
tion rates per unit area, with the median evolving from 
E* ~ 1.5 - 3.0 X IQ-'^Mq yr"^ kpc"^ at z = 4 - 5 to 
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Figure 11. The top left panel (a) shows the evolution of global volume weighted cosmic star density. Panel (b) shows the evolution of 
the global volume weighted star formation density. The lower panels (c) and (d) show the same quantities but weighted by cross-section 
(above the DL A threshold) . In each case we plot results for the three feedback models using the LC93 age distributions and for the coeval 
evolution model with the normal feed back model. The shaded regions in the figures show the volume averaged stellar densities and star 
formation rates from iHopkins I J2004) and iHopkins. Rao fc Turnshek I <2005D that encompass the majority of the observational points 
on the 'Madau'-diagram after cor rection for dust obscuration. The points plotted in panel (d) show the star formation rates inferred by 
iHopk ins. Rao &: Turnshek I JSoO^) for DLA systems (see text for details). 



our model, galaxies 

at high redshift are experiencing bursty star formation (with 
high star formation rates per unit area) followed by more 
quiescent star formation evolution at late times, when the 
gas infall rate is regulated by the Hubble time rather than 
the free-fall time (see EOO for details). The median values 
of S* ~ 10^'^ — W~^Mq yr~^ kpc ~^ are in good agree ment 
with the semi-analytic model of Okoshi et al. f2005') and 
with the numerical simulations of Naeaminc ct al. (2004t|). 

The evolution of the cosmological stellar density and 
star formation rates are shown in Fig. 1111 We show results 
for the three feedback models for the LC93 merger history, 
and for a model with coeval evolution. It is important to note 
that the volume averaged stellar density and star formation 
rates shown in this figure are weighted strongly towards the 



inner regions of massive disc systems, which make little con- 
tribution to the DLA cross-sections. As we will show below, 
these global volume averaged properties are largely discon- 
nected from the properties of the DLA population. It is also 
wor th noting that, unlike more elaborate se mi-analytic mod- 
els (|Croton et alJ 1200^1: ISower et al.l 120061 ). our model has 
not been fine-tuned to match these global properties. Nev- 
ertheless, it performs surprisingly well given its simplicity. 

All models show high star formation rates at high red- 
shift arising from the initial 'bursty' phase of star formation 
model, followed by a gentle decline in the comoving star for- 
mation density. The differences between the feedback mod- 
els are larger at higher redshifts, where the weaker feedback 
in the WFB model results in high star formation rates. At 
lower redshifts {z ^ 2) all three feedback models produce 
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Figure 12. Evolution of tlic cosmic hydrogen density in DLAs for our models ovcrplottcd with obs ervati ons from lZwaan et al 

(open square at z = 0), Rao ct al. (2006) (solid squares at z ~ 1) and Prochaska fc Hcrbcrt-For3 l2004l) (solid circles at z > 2). The 

model results have been corrected for a hydrogen fraction of Xh =0.7. 



star formation rates which are within a factor of two of each 
other. The global volume averaged star formation rates are 
quite sensitive to the assumed age distributions, especially 
at high redshift. 

Fig. Illb show the evolution of the volume aver- 
aged star formation rates in our models. The shaded re- 
gions in Figs 11a and lib shows t he range of obs erva- 
tional estimates from the papers of iHopk ins I teOO^Tl and 
[Hopkins, Rao & Turnshek I i2005() . including corrections for 
dust obscuration. Figs 11a and lib give an impression of the 
uncertainties in both the model predictions and the obser- 
vations. The models broadly match the extinction corrected 
star formation rates, though it is clear that a more accu- 
rate treatment of the merger history is required to produce 
reliable predictions. 

The main point of this comparison is illustrated in Figs. 
Illfc and llll i. These show the evolution of the cosmic stel- 
lar density and star formation rate but now weighted by 
cross-section above the DLA threshold, rather than volume 
averaged over all systems. According to our model, DLA 
systems account for a small fraction of the total star for- 
mation at all redshifts. This is broad l y con sistent with the 
results of lHopkins. R.ao fc Turnshek I (|200,4) . who infer the 
contribution of DLAs to fi* and dp»/dt by integrating over 
the observ e d DL A column-density distribution assuming a 
lKenniciit3 1^23) relation between t he star for mation rate 
and th e gas surface density. (The Hopkins. Rao fc Turns hekl 
l|200,'Th results are shown by the points in Fig. lid). Most 
of the star formation occurs in the inner parts of galaxies, 
whereas the DLA are sampling the outer gaseous discs where 



there is relatively little star formation. This, of course, is why 
the metallicity distributions are skewed to low metallicities 
at all redshifts. One should therefore be extremely cautious 
in interpreting models of cosmic chemical evolution that link 
the volume averaged star formation rates to the metallicity 
evolu t ion of DLA sel e cted systems dWolfe and Prochaska I 
119981: iPei et alJ Il999l : IWolfe et alJ |2005^ . As our models 
show, these quantities are linked only indirectly and are very 
likely probing star formation in very different environments. 



4.5 Cosmological evolution of Hhi 

From our models we can calculate the cosmological evolu- 
tion of neutral hydrogen. This is an important consistency 
check of the models, since unlike the stellar density, most 
of the volume averaged density in neutral hydrogen in our 
models is contained in DLA system s. Fig. 1121 s hows obser- 
vational results from |Z waan et al.ll2005a : .Rao et alJl2006l : 
IProchaska fc Herbert-Forn200^ spanning the entire red- 
shift range z ^ 5. The observations show little evidence 
for any evolution in ^hi over the redshift range 1 ^ 2 ^ 5. 

Our model predictions are also shown after correcting 
for a hydrogen fraction oi Xh = 0.7. No correction has been 
made for the fraction of hydrogen in molecular form. The 
results are relatively insensitive to the strength of feedback, 
but are quite sensitive to the assumed age distributions. The 
models with the LC93 age distributions underpredict Qhi 
at z ^ 2, whereas the coeval model agrees quite well with 
the observations over the redshift range 1 ^ z ^ 5. As with 
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several other properties discussed in this Section, the coeval 
models with those with the LC93 age distribution roughly 
bracket the observations. All of the models are high (by a 
factor of ~ 3) compared to the observed HIjiensity at z = 
(determined from 21 cm observations (iZwaan et al..,2 005a)'l. 
We do not regard this discrepancy as particularly serious, 
since it can be explained if some of the gas in the extended 
outer parts of small disc systems is stripped and returned to 
the intracluster or intra-group medium as structure grows. 

Most of the neutral hydrogen density in our models is in 
the outer parts of dwarf gala:xies and plays little role in the 
volume averaged star formation history of the Universe. The 
evolution of Qm at high redshift is therefore more closely 
linked to feedback processes and the spatial distribution of 
gas in the outer parts of dwarf galaxies than it is to the cos- 
mic star formation history. Our models show that the low 
metallicities, characteristic of DLAs, and the lack of evo- 
lution of Qhi are consistent with a picture in which the 
DLA sight-lines are preferentially sampling the outer parts 
of dwarf galaxies. 



5 CONCLUSIONS 

In this paper we have used a model of self-regulated 
star formation and a physically motivated model of stellar 
feedback to construct a simple model of damped Lyman-a 
systems. 

Our models reproduce the low mean metallicities seen 
in DLAs, the lack of evolution of the mean metallicity and 
can account, at least qualitatively, for the observed spread 
in metallicities at each redshift. In our model, most DLA 
sight-lines probe the outer gaseous parts of 'dwarf galaxies 
(i^vir ~ 70 kms~^), where the star formation rates and metal 
enrichment are always low. Occasionally, sight-lines intersect 
the inner regions of disc systems where the gas metallici- 
ties are high. Thus geometry is primarily responsible for the 
large spread in metallicities seen in the observations. The 
metallicity distributions in our models are relatively insen- 
sitive to the strength of the stellar feedback assumed and 
to the assumed merger histories of the galaxies (parame- 
terised by their age distributions). We therefore believe that 
the DLA metallicities are a robust feature of our model. 
The metallicities predicted by our model are lower, and in 
much better agreement wi th observations, than those found 
in num erical simulations iCen et al.l 120031 : iNaeamine et alJ 
|20043). Differences in the models for stellar feedback and 
star formation are the most plausible reasons for this. 

In many respects, our conclus ions are similar to those 
of previous semi-analyti c model s ("Okoshi ct al."20 04|) and 
numer ical simulations iCen et al. .2003, ; .Naeamine et alJ 
l2004d) . However, our model differs from previous semi- 
analytic calculations in that the infall model, star formation 
prescription and feedback model, all key factors in deter- 
mining the metallicities and spatial distributions of gas and 
stars, are all physically motivated. The model has, therefore, 
few adjustable parameters. Indeed, the model was developed 
specifically to investigate the role of stellar feedback during 
galaxy formation (EOO), yet with no further modification 
provides a reasonable description of the properties of DLA 
systems. 

Nevertheless, it is clear that the physics behind the DLA 



population is complex, and it is important to understand 
both the strengths and weaknesses of a simple semi-analytic 
model of the type described here. Although the metallicity 
predictions are robust, the sizes of the gas discs in our model, 
particularly at high redshifts, are sensitive to the assumed 
age distributions. These distributions, in turn, depend on 
the merger history of the gaseous discs. To illustrate the 
sensitivity of our results to the age distributions, we have 
used two simple models, one based on the age distributions 
of the dark haloes (using the prescription of LC93) and the 
other based on the assumption that all systems formed at 
t = (the coeval model). The actual age distribution for 
the baryonic components is likely to lie between these two 
extremes. 

The rates-of-incidence (Figure 8) are particularly sensi- 
tive to the total gas cross-section, and hence to the model of 
the age distribution. Our model cannot predict this quan- 
tity reliably. Nevertheless, the two models for the age distri- 
butions roughly bracket the observations over the redshift 
range 1.5 z ^ 5, suggesting that a more realistic merger 
model may be able to account for the data over this redshift 
range. In fact, the numerical simulations of 'Nag amine et al] 
^2004a^ reproduce the observed rates-of-incidence quite well, 
and (not surprisingly) their cross-sections are intermedi- 
ate between those of our two age-distribution models. (The 
physical mechanisms governing the cross-sections in the nu- 
merical simulations are not clear, however.) At lower red- 
shift, both age-distribution models overpredict the rates- 
of-incidence, suggesting that we are missing some physical 
mechanism that can limit the growth of extended gaseous 
discs at low redshift. The cosmic density in neutral hydro- 
gen is also dependent, but to a much lesser degree, on the 
assumed age distributions. As with the rates-of-incidence, 
the two models for the age distribution roughly bracket the 
observations, except at z = 0, where our models overpredict 
the local HI density by a factor of ~ 1.5-2.5. 

The relation between the metallicities and column den- 
sities is not well reproduced by our models, particularly at 
redshifts z ^ 2. In our models, the gaseous discs at any time 
have roughly constant surface density and truncate abruptly 
(with an outer radius fixed by the angular momentum of in- 
falling gas). As a result, our models fail to reproduce the 
spread to low column densities seen in the observations 
(though only by a factor of two or so) . The abrupt trunca- 
tion is clearly an artificial feature of our models and it is easy 
to think of a number of physical mechanisms (discussed in 
Section 3.2) that would lead to lower column densities in the 
outer parts of galaxies, without altering the metallici t y of the 
gas. The numerical simulations of Nagami ne et al.l i2004lj) 
do show a large range of column densities, extending below 
the DLA threshold. The observed high proportion of low 
metallicity, high column density density (A^'hi > 10^^ cm~^) 
systems at z 2 is more problematic, both for our model 
and the numerical simulations. The resolution of this dis- 
crepancy is not at all clear. 

Our models can roughly reproduce the observed volume 
averaged stellar density, HI density and star formation rate 
as a function of redshift (though the models were not finely 
tuned to do so). The cosmic density in HI is dominated 
by gas in DLA systems. However, the stellar density and 
star formation rates weighted by cross-section are consider- 
ably lower than their cosmic values over the entire redshift 
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range ^ 2 5. In our models, the galaxies responsible for 
DLAs make only a small contribution to the total stellar den- 
sity produced in the Un iverse. This conclusion is broadly i n 
line with th e resu lts of Hopkins. Ra o fc Turnshek I i2005l) . 
who used a iKennicutt. 1.1998,') relation to infer the contri- 
bution of DLA systems to the volume averaged star forma- 
tion rate. According to our model, there is only a indirect 
relation between the volume averaged star formation rate, 
flm and DLA metallicities. Thus one should be skeptical 
of models of 'cosmic' chemi cal evolution that attempt t o 
relate these qu a ntities (e.g. IWolfe and Prochaska I il99d) : 
iPei et alJ lUggl l : I Wolfe et alj tOO^ ). 

As we have mentioned in the Introduction, one of the 
main purposes of a semi-analytic model is to gain physi- 
cal insight into complex problems. Ultimately, given high 
enough spatial resolution, it should be possible to develop 
realistic numerical hydrodynamic models of DLAs. In the 
meantime, we pose the following problems for numerical sim- 
ulators: 

(i) In our model, the sizes (and hence cross-sections) of 
the DLA systems are set by the high angular momentum 
gas within the dark haloes, assuming angular momentum is 
strictly conserved. (The low angular momentum gas is pref- 
erentially expelled in a wind). Is this really true? 

(ii) We have found that some properties of the DLA pop- 
ulation are sensitive to the age distribution and hence the 
merger history of the baryonic systems within dark haloes. 
What is the actual merger history of these baryonic sys- 
tems? To what extent does this merger history affect the 
properties of DLAs. 

(iii) How sensitive are the metallicities of DLAs to star for- 
mation and stellar feedback? Is it possible to incorporate 
more realistic models of these processes in numerical simu- 
lations? Is it possible to resolve the discrepanci es between 
observations and simulations performed so far JCen et alJ 
l2003l:lNagamine et aljEo04bl) ? 

(iv) Is is possible to reproduce the observed proportion of 
low metallicity systems with high column density, that nei- 
ther our model or the numerical simulations seem able to 
explain? 

(v) Clearly a theoretical model in which dark haloes are 
populated by single circular, homogeneous, gaseous discs is 
a gross over-simplification. To what extent are tidal features 
and other irregularities important in understanding the DLA 
population? How important are inhomogeneities, such as in- 
dividual star forming regions, in explaining the high end of 
the metallicity distribution? 

(vi) Our model is too simplified to model the velocity struc- 
ture seen in metal-line-systems associated with DLAs. Can 
this structure be reproduced by detailed hydrodynamical 
simula tions, as suggested by the work by iHaehnelt et alJ 
Jl998h ? 

(vii) In our model, the DLAs make a sub-dominant contri- 
bution to the cosmic star formation rate at all redshifts. Is 
this really true? 

No doubt the reader can think of many other prob- 
lems. Despite the compelling case for further work, we do 
beheve that the simple model presented here gives a plausi- 
ble explanation for the observed metallicities of DLAs, and 
that these metallicities can only be understood if the DLAs 



are predominantly 'dwarf-like systems that contribute lit- 
tle to the net cosm ic stellar density. The 'classical' picture 
JWolfe et alJlTgS^ of DLAs as giant spiral discs slowly con- 
verting their gas into stars seems to us to be irreconcilable 
with the observations. Furthermore, we now know enough 
about the primordial cosmological fluc tuations, from Mpc 
scales to the scale of the Hubble radius JSoereel et alJl2003l : 
ICole et al.l2'005tlViel et al. l2004ls'eliak et aL l2005l) that it 
is difficult to imagine a cosmogonic context for the classical 
picture. 
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